
clear all
set more off


cap cd "E:\Dropbox\paper_returnmig\estimation\replication_files/"
************************************* Stylized facts
* Figure 1 requires downloading the data from the Migration Data Portal: https://www.migrationdataportal.org/
/*
twoway  (connect pending year, lc(black) lp(dash) mc(black) msym(Oh)) ///
		(connect returnees year, yaxis(2)  lc(gs6) lp(solid) mc(gs6) msym(Dh) mc(gs10) ///
		ytitle("Asylum seekers (in 1,000)", axis(1)) ytitle("Assisted returnees (in 1,000)", axis(2)) ///
		ylabel(, nogrid) ///
		legend(order(1 "Pending asylum abroad" 2 "AVR & resettlement returnees") rows(2)))
graph export "figure_1.pdf", replace
*/
*************************************


cap use "Auer_Schaub_senegambia.dta", clear
*************************************
global red "200 0 0"
global blue "68 119 170"
global yellow "204 187 68"
global green "34 136 51"
global orange "233 118 0"
global magenta "150 79 142"
**************************************
global mainsample "if region != 2"
global dakarsample "if region == 2"
global indiv female age single i.education  
global covars job ib(2).hhecon hhland resp_ret_eur resp_ret_afr   fam_mig_any  fam_ret_any
global fe pop_all i.location  
global localmig mig_any
global covars_lifesat	female age single i.education ib(2).lifesat resp_ret_eur resp_ret_afr   fam_mig_any  fam_ret_any
global covar_econ 		female age single i.education hhland resp_ret_eur resp_ret_afr
**************************************



************************
* descriptives on reason for return
g ret_involuntary = .
	replace ret_involuntary = 0 if resp_ret_eur == 1 & retres == 0
	replace ret_involuntary = 1 if resp_ret_eur == 1 & retres != 0
tab ret_involuntary
tab retres if resp_ret_eur == 1


********************************************************************************
* Summary Statistics
reg emigplan_eur ret_eur $localmig $indiv $covars $fe $mainsample, robust 
g insample=e(sample)

eststo:  estpost summarize female age single educ1 educ2 educ3 educ4 job hhecon hhland resp_ret_eur resp_ret_afr ///
							pop_all location_female location_islam location_young district_famsize ///
							fam_mig_any fam_ret_any mig_afr mig_eur ret_afr ret_eur ///
							if insample == 1
esttab using "table_a1.tex", append cells("mean(fmt(3)) sd(fmt(3)) min(fmt(3)) max(fmt(3))") label nodepvar 


********** endogenous return economic factors?
*********************** Variation in migration perception by whether returnees are present
foreach i in female age single educ1 educ2 educ3 educ4 job hhecon hhland resp_ret_eur resp_ret_afr fam_mig_any  fam_ret_any pop_all    {
	reg bin_ret_any `i'  $localmig  i.location, robust  cluster(location)  
	estimates store `i'
	reg bin_ret_any `i' , robust cluster(location)
	estimates store `i'_1	
}

global color2 gs10 gs10
global offset 0.35
coefplot (female, msym(Oh) ciopts(lcolor(black))) (female_1, msym(Dh) mc(gs10) ciopts(lcolor($color2)) offset($offset)) ///
		(age, msym(Oh) ciopts(lcolor(black))) (age_1, msym(Dh) mc(gs10) ciopts(lcolor($color2)) offset($offset)) ///
		(single, msym(Oh) ciopts(lcolor(black))) (single_1, msym(Dh) mc(gs10) ciopts(lcolor($color2)) offset($offset)) ///
		(educ1, msym(Oh) ciopts(lcolor(black))) (educ1_1, msym(Dh) mc(gs10) ciopts(lcolor($color2)) offset($offset)) /// 
		(educ2, msym(Oh) ciopts(lcolor(black))) (educ2_1, msym(Dh) mc(gs10) ciopts(lcolor($color2)) offset($offset)) ///
		(educ3, msym(Oh) ciopts(lcolor(black))) (educ3_1, msym(Dh) mc(gs10) ciopts(lcolor($color2)) offset($offset)) ///
		(educ4, msym(Oh) ciopts(lcolor(black))) (educ4_1, msym(Dh) mc(gs10) ciopts(lcolor($color2)) offset($offset)) ///
		(job, msym(Oh) ciopts(lcolor(black))) (job_1, msym(Dh) mc(gs10) ciopts(lcolor($color2)) offset($offset)) ///
		(hhecon, msym(Oh) ciopts(lcolor(black))) (hhecon_1, msym(Dh) mc(gs10) ciopts(lcolor($color2)) offset($offset)) ///
		(hhland, msym(Oh) ciopts(lcolor(black))) (hhland_1, msym(Dh) mc(gs10) ciopts(lcolor($color2)) offset($offset)) ///
		(fam_mig_any, msym(Oh) ciopts(lcolor(black))) (fam_mig_any_1, msym(Dh) mc(gs10) ciopts(lcolor($color2)) offset($offset)) ///
		(fam_ret_any, msym(Oh) ciopts(lcolor(black))) (fam_ret_any_1, msym(Dh) mc(gs10) ciopts(lcolor($color2)) offset($offset)) ///
		(pop_all, msym(Oh) ciopts(lcolor(black))) (pop_all_1, msym(Dh) mc(gs10) ciopts(lcolor($color2)) offset($offset)) ///
		, legend(order(6 "Before adjustment" 3 "After adjusting for region, current migrants") rows(2) size(small)) xline(0, lc(red)) nolabel drop(_cons  mig_afr mig_eur ?.location) ///
	 msym(Oh Dh) xlabel(-0.75(0.25)0.75, labsize(small)) levels(95 90) grid(none) /// 
	 ylabel(,labsize(small)) ///
	rename(female = "Female" female_1 = "Female" age = "Age" age_1 = "Age" single = "Single" single_1 "Single" ///
			educ1 = "No education" educ1_1 "No education" educ2 = "Primary education" educ2_1 "Primary education" ///
			educ3 = "Lower secondary" educ3_1 "Lower secondary" educ4 = "Upper secondary +" educ4_1 "Upper secondary +" ///
			job = "Paid work" job_1 "Paid work" hhecon = "HH economic situation" hhecon_1 "HH economic situation" ///
			hhland = "HH land ownership" hhland_1 = "HH land ownership" ///
			pop_all = "Local population" pop_all_1 = "Local population" ///
			fam_mig_any = "Family migrants" fam_mig_any_1 = "Family migrants" ///
			fam_ret_any = "Family returnees" fam_ret_any_1 = "Family returnees" ///
			location_1 = "Gambia" location_1_1 "Gambia" ///
			location_2 = "Kolda" location_2_1 "Kolda" ///
			location_3 = "Sedhiou" location_3_1 "Sedhiou" ///
			location_4 = "Ziguinchor" location_4_1 = "Ziguinchor")
	graph export "figure_5.pdf", replace
	
global distrespvars av_age av_female av_education av_job av_hhecon

foreach i in  $distrespvars  {
	reg anyret_dist `i'  mig_eur_dist  mig_afr_dist if distid == 1  ,  robust cluster(location)
	estimates store `i'
	reg anyret_dist `i'   if distid == 1 & mig_any_dist != . , robust cluster(location)
	estimates store `i'_1	
}
	
foreach i in  location_female  location_young location_islam {
reg anyret_dist `i' mig_afr_dist mig_eur_dist  if distid == 1 
	estimates store `i'
reg anyret_dist `i'   if distid == 1 & mig_any_dist != .
	estimates store `i'_1
}

global offset 0.15
coefplot (location_female, msym(Oh) ciopts(lcolor(black))) (location_female_1, msym(Dh) mc(gs10) ciopts(lcolor($color2)) offset($offset)) ///
		(location_young, msym(Oh) ciopts(lcolor(black))) (location_young_1, msym(Dh) mc(gs10) ciopts(lcolor($color2)) offset($offset)) ///
		(location_islam, msym(Oh) ciopts(lcolor(black))) (location_islam_1, msym(Dh) mc(gs10) ciopts(lcolor($color2)) offset($offset)) ///		
		(av_age, msym(Oh) ciopts(lcolor(black))) (av_age_1, msym(Dh) mc(gs10) ciopts(lcolor($color2)) offset($offset)) ///
		(av_female, msym(Oh) ciopts(lcolor(black))) (av_female_1, msym(Dh) mc(gs10) ciopts(lcolor($color2)) offset($offset)) ///
		(av_education, msym(Oh) ciopts(lcolor(black))) (av_education_1, msym(Dh) mc(gs10) ciopts(lcolor($color2)) offset($offset)) ///
		(av_job, msym(Oh) ciopts(lcolor(black))) (av_job_1, msym(Dh) mc(gs10) ciopts(lcolor($color2)) offset($offset)) ///
		(av_hhecon, msym(Oh) ciopts(lcolor(black))) (av_hhecon_1, msym(Dh) mc(gs10) ciopts(lcolor($color2)) offset($offset)) ///
	, yline(3.5, lc(black) lp(solid) lw(vthin)) legend(order(6 "Before adjustment" 3 "After adjusting for current migrants") rows(2) size(small)) xline(0, lc(red)) nolabel drop(_cons mig_eur_dist mig_afr_dist ?.location) ///
	 msym(Oh Dh) xlabel(, labsize(small)) levels(95 90) grid(none) /// 
	 ylabel(,labsize(small)) xlabel(,labsize(small)) /// -15(5)15
	 rename(location_female = "Female share in district" location_female_1 = "Female share in district" ///
	 location_young = "Young residents share in district" location_young_1 = "Young residents share in district" ///
	 location_islam = "Muslim share in district" location_islam_1 = "Muslim share in district" ///
	 av_age = "Mean respondent age" av_age_1 = "Mean respondent age" ///
	 av_female = "Share female respondents" av_female_1 = "Mean female respondents" ///
	 av_education ="Mean respondent education" av_education_1 = "Mean respondent education" ///
	 av_job = "Share employed respondents" av_job_1 = "Share employed respondents" ///
	 av_hhecon = "Mean HH economic situation" av_hhecon_1 = "Mean HH economic situation")
graph export "figure_4.pdf", replace


	
********************************************************************************
********** MAIN ANALYSIS: CONDITIONAL ON MIGRANTS and POP, RETURNEES outside family are AT RANDOM
reg emigplan_afr ret_afr $localmig $indiv $fe $mainsample, robust 
	estimates store m1

reg emigplan_afr ret_afr $localmig $indiv $covars $fe $mainsample, robust
	estimates store m2	
	
reg emigplan_afr ret_afr ret_eur $localmig $indiv $covars $fe $mainsample, robust 
	estimates store m3
	
reg emigplan_eur ret_eur  $localmig $indiv $fe $mainsample, robust
	estimates store m4

reg emigplan_eur ret_eur  $localmig $indiv $covars $fe $mainsample, robust 
	estimates store m5

reg emigplan_eur ret_eur ret_afr $localmig $indiv $covars $fe $mainsample, robust 
	estimates store m6
	
esttab m1 m2 m3 m4 m5 m6 using "table_1.tex", replace label b(3) se(3) stats(r2 N) par plain star(* 0.10 ** 0.05 *** 0.01)
est clear



**************** MECHANISM I: Material consequences
***** A. HH Economic situation
reg econ_hh z1_ret_afr $covar_econ z1_mig_afr z1_mig_eur  $fe   $mainsample, robust 
	estimates store t1

reg econ_hh z1_ret_afr $covar_econ z1_mig_afr z1_mig_eur z1_fam_mig_afr z1_fam_ret_afr $fe   $mainsample, robust 
	estimates store t2
	
	
reg econ_hh z1_ret_eur $covar_econ z1_mig_eur z1_mig_afr $fe   $mainsample, robust 
	estimates store t3
	
reg econ_hh z1_ret_eur $covar_econ z1_mig_eur z1_mig_afr z1_fam_mig_eur z1_fam_ret_eur $fe   $mainsample, robust 
	estimates store t4	

	
esttab t1 t2 t3 t4 using "table_2.tex", replace label b(3) se(3) stats(r2 N) par plain star(* 0.10 ** 0.05 *** 0.01) keep(z1_ret_* z1_mig_* z1_fam_mig_* z1_fam_ret_*)
est clear

***** B. HH remittances as main income
g remittances = 0
	replace remittances = 1 if hhsector == 4

reg remittances z1_ret_afr $covar_econ z1_mig_afr z1_mig_eur  $fe   $mainsample, robust 
	estimates store t1

reg remittances z1_ret_afr $covar_econ z1_mig_afr z1_mig_eur z1_fam_mig_afr z1_fam_ret_afr $fe   $mainsample, robust 
	estimates store t2
	
	
reg remittances z1_ret_eur $covar_econ z1_mig_eur z1_mig_afr $fe   $mainsample, robust 
	estimates store t3
	
reg remittances z1_ret_eur $covar_econ z1_mig_eur z1_mig_afr z1_fam_mig_eur z1_fam_ret_eur $fe   $mainsample, robust 
	estimates store t4	

	
esttab t1 t2 t3 t4 using "table_2.tex", append label b(3) se(3) stats(r2 N) par plain star(* 0.10 ** 0.05 *** 0.01) keep(z1_ret_* z1_mig_* z1_fam_mig_* z1_fam_ret_*)
est clear	
	
	
	
**************** MECHANISM II: Perceptions

global percepvars neglife_bin ret_fail_bin
su $percepvars
global localmigrants z1_mig_afr z1_mig_eur 

foreach i in $percepvars {	
if "`i'" == "neglife_bin"	{
	global storing replace
}
else {
	global storing append
}	
reg  `i'  z1_ret_afr $indiv job ib(2).hhecon hhland resp_ret_eur resp_ret_afr $localmigrants $fe   $mainsample, robust 
	estimates store t1

reg  `i'  z1_ret_afr $indiv job ib(2).hhecon hhland resp_ret_eur resp_ret_afr  $localmigrants z1_fam_mig_afr z1_fam_ret_afr $fe   $mainsample, robust
	estimates store t2
	
reg  `i'  z1_ret_eur $indiv job ib(2).hhecon hhland resp_ret_eur resp_ret_afr $localmigrants $fe   $mainsample, robust
	estimates store t3
	
reg  `i'  z1_ret_eur $indiv job ib(2).hhecon hhland resp_ret_eur resp_ret_afr $localmigrants z1_fam_mig_eur z1_fam_ret_eur $fe   $mainsample, robust
	estimates store t4	
	
esttab t1 t2 t3 t4 using "Table_3.tex", $storing label b(3) se(3) stats(r2 N) par plain star(* 0.10 ** 0.05 *** 0.01) keep(z1_ret_* z1_fam_mig_* z1_fam_ret_*)
est clear
}

************** ROBUSTNESS: show that populations in districts are homogeneous. we adjust for pop in district to account to minor variations.
forvalues i = 0/4 {
   if `i' == 0 {
       global title "The Gambia (GMB)"
   }
   else if `i' == 1 {
       global title "Kolda (SEN)"
   }
   else if `i' == 2 {
       global title "Sedhiou (SEN)"
   }
   else if `i' == 3 {
      global title "Ziguinchor (SEN)"
   }   
   kdensity pop_all if location == `i', bw(100) name(g`i', replace) nor  title("$title", size(medsmall)) note("") ///
		xlabel(0(500)2000, nogrid) xtitle(" " "Local population") ylabel(0(0.001)0.002, nogrid) ///
		legend(order(1 "Observed regional distribution" 2 "Normal distribution") rows(1) size(small))
}
grc1leg  g0 g1 g2 g3
graph export "figure_a1.pdf", replace	
	
*********************** Coorelation between current migrants and reutrnees
global size huge

twoway (scatter ret_any mig_any $mainsample & distid == 1, msym(oh) msize(large) mc(black%75) jitter(7) name(g1, replace) ///
		legend(off) xtitle("A. Current migrants from area", size($size)) ytitle("Local returnees", size($size)) xlabel(0(20)60, labsize($size)) ylabel(0(5)15, nogrid labsize($size)))

twoway 	(scatter av_plans mig_any $mainsample & distid == 1, msym(oh) msize(large) mc(black%75) jitter(7)  name(g2, replace) ///
		legend(off) xtitle("B. Current migrants from area", size($size)) ytitle("Share with migration plans", size($size)) xlabel(0(20)60, labsize($size)) ylabel(0(25)100, nogrid labsize($size))) 

twoway 	(scatter av_plans ret_any $mainsample & distid == 1, msym(oh) msize(large) mc(black%75) jitter(7) name(g3, replace) ///
		legend(off) xtitle("C. Local returnees", size($size)) ytitle("Share with migration plans", size($size)) xlabel(0(5)15, labsize($size)) ylabel(0(25)100, nogrid labsize($size)))

graph combine g1 g2 g3, rows(1) xsize(7) ysize(2)
graph export "figure_a2.pdf", replace


********************************************************************************
************** ROBUSTNESS: 
* 1: replace economic variables with life satisfaction variables in main model  
reg emigplan_afr ret_afr  $localmig  $covars_lifesat $fe $mainsample, robust  
	estimates store m1

reg emigplan_eur ret_eur $localmig  $covars_lifesat $fe $mainsample, robust 
	estimates store m2
	

* 2: use standardized IVs	
egen std_ret_eur = std(ret_eur) $mainsample & emigplan_any != . 
egen std_ret_afr = std(ret_afr) $mainsample & emigplan_any != . 

reg emigplan_afr  std_ret_afr $localmig $indiv $covars $fe $mainsample, robust 
	estimates store m3

reg emigplan_eur std_ret_eur  $localmig $indiv $covars $fe $mainsample, robust 
	estimates store m4
	
esttab m1 m2 m3 m4 using "table_a2.tex", replace label b(3) se(3) stats(r2 N) par plain star(* 0.10 ** 0.05 *** 0.01)
est clear


****** FE by geographical location
forvalues i = 1/5 {
eststo a`i': reg emigplan_afr ret_afr  $indiv $covars i.geo_`i' $mainsample, robust 
eststo b`i': reg emigplan_eur ret_eur  $indiv $covars i.geo_`i' $mainsample, robust
}

eststo a6: reg emigplan_afr ret_afr  $indiv $covars i.district_num $mainsample, robust
eststo b6: reg emigplan_eur ret_eur  $indiv $covars i.district_num $mainsample, robust


esttab a1 a2 a3 a4 a5 a6 using "table_a3.tex", replace label b(3) se(3) stats(r2 N) par plain star(* 0.10 ** 0.05 *** 0.01) keep(ret_afr)
esttab b1 b2 b3 b4 b5 b6 using "table_a3.tex", append label b(3) se(3) stats(r2 N) par plain star(* 0.10 ** 0.05 *** 0.01) keep(ret_eur)
distinct geo_1 geo_2 geo_3 geo_4 geo_5 district_num $mainsample & female != . & insample == 1
est clear

	
	
********** OUT fo Sample test
reg emigplan_afr ret_afr $localmig $indiv $fe $dakarsample, robust 
	estimates store m1

reg emigplan_afr ret_afr $localmig $indiv $covars $fe $dakarsample, robust  
	estimates store m2

	
reg emigplan_eur ret_eur $localmig $indiv $fe $dakarsample, robust 
	estimates store m3

reg emigplan_eur ret_eur $localmig $indiv $covars $fe $dakarsample, robust 
	estimates store m4
	
esttab m1 m2 m3 m4  using "table_a4.tex", replace label b(3) se(3) stats(r2 N) par plain star(* 0.10 ** 0.05 *** 0.01)
est clear

	
********** OMITTED VARIABLE BIAS
* OSTER ROBUSTNESS TEST
reg emigplan_eur ret_eur  $mainsample, robust
reg emigplan_eur ret_eur $localmig $indiv $covars $fe $mainsample, robust 
psacalc delta ret_eur, rmax(0.05174)
psacalc delta ret_eur, rmax(0.25)

* SENSEMAKR
reg emigplan_eur ret_eur $localmig $indiv $covars $fe $mainsample, robust 
reg emigplan_eur ret_eur $localmig $indiv $covars $fe $mainsample,  
reg emigplan_eur ret_eur   

sensemakr emigplan_eur ret_eur $localmig $indiv $covars $fe $mainsample,  ///
         treat(ret_eur) benchmark(job) kd(10 20) ky(10 20)  latex(sensmakr)  contourplot
graph export "figure_a3.pdf", replace	
		 